# TODO: Add comment
# 
# Author: guochun
###############################################################################


#given r

obf=function(x,r,k,q){
	#x=c("lambda","a")
	if(any(x<0))
		return(-Inf)
	lam=x[1]
	a=x[2]
	n=length(r)
	#re=n*log(lam)+sum(log(r))-(k+1)*sum(log(1+(pi*lam*r^2)/k))
	re=n*log(lam)+n*lgamma(k+a)+n*a*log(a)-n*lgamma(a)+
			sum(log(r)+(k-1)*log(pi*lam/q*r^2)-(k+a)*log(pi*lam/q*r^2+a))
	return(re)
}

pSizeTor.neg=function(samples){
	sampleMethod=attr(samples,"sampleMethod")
	spclass=class(sampleMethod)
	if(spclass=="angleOrder"){
		q=sampleMethod@q
		k=sampleMethod@k
		if(sampleMethod@composite){
			r=unlist(samples)
		}else{
			r=as.numeric(samples)
		}
		del=which(is.na(r) | r==Inf)
		if(length(del)!=0)
			r=r[-del]
		if(length(r)>0){
			ngpars=optim(c(0.3,5),fn=obf,r=r,k=k,q=q,control=list(maxit=100000,fnscale=-1))
			if(ngpars$convergence==0)
				N=ngpars$par[1]*500*1000
			else
				N=NA
		}else{
			N=NA
		}
		
		return(N)
	}
	return(NA)
}






